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NOTATION 


The following symbols art used in this paper : 

A ■ cyclic strength coefficient; 

£ - expanding point of the Taylor's series; 
b ■ fatigue strength exponent; 
c - fatigue ductility exponent; 

E « modulus of elasticity; 
f - probability density function; 
g(Jg)« failure function of Jg; 
h(*)“ function defined by Eq. 21; 
k^- stress concentration factor; 
k » number of normal equations; 
m - nunber of observations or data sets; 

N, N ■ cycles to failure , service life in cycles; 
n * number of random design variables; 
p(*)“ probability of; 

P£ a probability of failure; 

S ■ nominal stress; 
s ■ cyclic strain hardening exponent; 

]g ■ design parameter vector; 
uj« defined by Eq. 29; 

Uj" reduced variable; 
desing point values; 

A 

a, a ■ coefficient in polynomial equation, least square estimator 
0 - safety index; 


i 



c, ■ strain rang*, cyclic strain amplitude; 

e| ■ fatigue ductility coefficient; 

X - direction cosine; 

o, o q - stress range, mean stress; 

'• fatigue strength coefficient; 

*(•)■ standard normal cdf; 
ft ■ failure domain. 
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EXECUTIVE SUMMARY 

Addressed herein Is the problem of computing the probability of failure 
for a structural component whose design factors are modelled as random vari- 
ables. The relationship between the design factors for a particular mode 
of failure is complicated; it can be defined only using a computer algorithm, 
e.g., finite element analysis, local strain analysis. The situation described 
relates especially to fatigue risk assessment in which a local strain analysis 
sucn as strain range partitioning'' is used to describe fatigue behavior of a 
component. But the procedure is general and can be applied to any mode of 
structural or mechanical failure. 

As an example, consider the fatigue problem described in the figure. 

A harmonic load (rero-to-peak) is applied to the part. The stress in the 
bulk of the material remains elastic, but there is cyclic plasticity at 
the notch. Local strain analysis, requiring a computer solution, can be 
employed to estimate cycles to fatigue crack initiation at the notch. Ran- 
dom design factors include (a) the stress amplitude S reflecting statistical 
scatter in load data and modelling errors, and (b) the fatigue strength co- 
efficient and fatigue ductility coefficient representing scatter in 
material fatigue behavior. Other sources of uncertainty could include, in 
general, modelling error in using Neuber's and Miner's Rule in the analysis 
and uncertainties in the cyclic stress strain curve. 

Cycles to crack initiation, N, is uncertain (a random variable) because 
it depends on other factors which are uncertain (random variables). The 

basic goal of a reliability analysis, . . . given the service life, N , 

o 

what is the probability of failure, i.e., the probability that N will be 

less than N ? 

o 
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A FATIGUE PROBLEM WHERE LOCAL STRAIN ANALYSIS CAN BE EMPLOYED TO PREDICT 

FATIGUE CRACK INITIATION 


Describes scatter in loading 



CYCLIC STRAIN-LIFE CURVE 


CYCLIC STRESS-STRAIN CURVE 



Stress 

Amplitude 

0 

a 
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Computer coats for a Monte Carlo solution would be excessive because 
the local strain program would have to be run for each random sample 
point (maybe 10,000 times). An alternative proposed herein suggests a 
simplified scheme for establishing an approximate algebraic relationship 
between the design variables. Then an advanced computational method 
(called the Rackwitc-Flessler algorithm) is used to estimate the proba- 
bility of failure. 

The procedure is much faster than Monte Carlo and probably more ac- 
curate for typical reliability problems although the latter has yet to 
be proven. In addition to fatigue reliability analysis, the method is 
thought to be particularly useful in any mechanical or structural problem 
where computer analysis is required to relate the random design variables. 
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Numerous methods have been developed during the past decade to de- 
scribe structural reliability without having to perform a multiple integration. 
For example, Hasofer and Lind [7] introduced the concept of generalized 
safety index, 8, as a "measure*' of reliability which used only mean and 
standard deviation of each design variable. However, unless the limit 
state is linear and all the random variables involved are normally distrib- 
uted, no direct relationship between reliability and the safety index 
can be described in the Hasofer-Lind formulation. Distributional informa- 
tion is not introduced. 

To extend the concept of the safety index, Rackwitz and Fiessler 
[12, 13] suggest that non-normal distributed variables can be transformed 
to equivalent normal variables in the Hasofer-Lind analysis. The refined 
safety index, 8, and the design points are obtained through a convergent 
numerical process. An approximation to the actual probability of 
failure is provided by letting p^ = $(-8), where 4> is the standard normal 
distribution function. A condition is that the failure surface is approxi- 
mately linear in the neighborhood of the design point. The theory behind 
this well-known R-F algorithm was later formally investigated by Ditlevsen 
[3] who called it the "principle of normal tail approximation." Chen and 
Lind [1] have recently refined the Rackwitz-Fiessler method by employing 
an additional scale parameter for the equivalent normal distribution. A 
more accurate approximation of the original distribution is promised. 

The performance of the R-F algorithm and the C-L algorithm, both 
methods of "fast probability integration," have been investigated [16]. 
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Both are considered to be effective and accurate for general reliability 
purposes in design. The method is fast in that numerical integration is 
avoided; convergence to the design point is rapid. However, the schema 
becomes very inefficient or very difficult to apply in some cases where 
the limit state cannot be expressed in a closed form equation. An alterna- 
tive is to use Monte Carlo simulation which is able to provide approximate 
solutions to any difficult problem. Unfortunately, Monte Carlo is costly, 
particularly when additional computation must be made to define the limit 
state. To compound the problem, very large sample sizes are required in 

order to get sufficient resolution in the tall areas to approximate the 

-3 -5 

low probabilities typical of structural application (e.g., 10 to 10 ). 

To overcome the difficulties associated with Monte Carlo and a full 
distributional approach, a numerical method is proposed whereby the limit 
state is approximated as a polynomial in the neighborhood of the design 
point. Then the fast probability integration technique can be applied 
effectively to estimate p^. 

DEFINITIONS AND PROBLEM FORMULATION 

Let JJF ■ (U^, Uj, • • • U^) denote the vector of design factors which 
are, in general, random variables. It will be assumed that all U^ are un- 
correlated with known distributions. The failure function , Z ■ g(JJ) of 
random design factors is formulated so that the event of failure is. 


The limit state is defined as. 
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gty) * 0 (2) 

The probability of failure is defined as, 

P f - P (g(JJ) 1 0] (3) 

P, can be evaluated from the following multiple integral 



where f (J^J) is the joint probability density function of the random vector, 

K 

II; ft is the failure domain in design parameter space. 

In general, evaluation of Eq. 4 requires numerical Integration. 

The operation is particularly difficult if the number of the random variables 
exceeds two. And in practice, the failure function does not necessarily 
have to be expressed in explicit form. An example presented herein is 
a computation of reliability for a fatigue problem for which the local 
strain approach is used to determine cycles to failure in a notched member to 
which a constant amplitude oscillatory load is applied. Fatigue craclc lnltla 
tion life N at the notch (point of stress concentration) which experiences 
cyclic plasticity is predicted using the following three equations [2, 6] 

1. Neuber's rule: 

(kS) 2 

oe - — (5) 

2. Cyclic stress-strain curve: 
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3. Strain-life curve which deli nee fatigue etrength. 


6 a " 


oi - a . o' - o c ^ b , 

<2N) b + ( f g r" £ ) ® f <2N) C 


where 

- streee concentration factor 

S - nominal atreae 

E » modulus of elasticity 

o - local stress 

c - local strain 

A - cyclic strength coefficient 

a - cyclic strain hardening exponent 

e # ■ cyclic strain amplitude 

o - mean stress 
o 

c* ■ fatigue strength coefficient 
e * » fatigue ductility coefficient 
b - fatigue strength exponent 
c • fatigue ductility exponent 


( 7 ) 


The details of computing fatigue initiation life is described by Collins 
[2] and by Fuchs and Stephens [6], The algorithm is complicated, and a 
computer program is required. Therefore, an algebraic expression of the 
failure function g(JJ) is not immediately available. 

To simplify the reliability problem, assume that only r£, end 
S (applied nominal stress range of a constant amplitude process) are random 
variables. The remaining parameters are constant. It will be assumed 
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also that the mean stress effect is negligible due to stress relaxation; 

o ■ 0. It is important to note that in the example which follows. It is 
o 

possible, with no extent ion of the theory, to Include, for example, 

, A, s, and even o q as random variables. 

A local strain analysis computer program is thfc only practical way of 
computing cycles to failure N as a function of , and S. A closed 

form expression of 

N - H (S , o’., e^) (8) 

does not exist. Given that the design life is N q , the event of failure 
is (N 1 N q ) and the limit state is (N ■ N q ) . 

Gqs. 5 to 7 could be expressed Implicitly as: 

S - h(o’ f , c' f , N), etc. (9) 

If the closed form expression of h were given the R-F algorithm could 
be used effectively to estimate given (N - N q ). The scheme presented 
here is to construct an approximation of h in the neighborhood of the de- 
sign point and then use the R-F method. Monte Carlo is used as a check 
on the R-F computations. 

MONTE CARLO METHOD 

To estimate probabilities associated with a complicated function of 
random variables, a Monte Carlo solution can be used. Consider 

Z " g( $ (10) 

where JJ is a vector of random variables. The computer generates a random 
sample* of ; i - 1, J, where J is typically a large number, e.g., 10,000. 
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Eq. XO Is employed to obtain a random sample of Z^ of slsa J. To estimate 
P(Z i 0), for example, one can take the number of Z^ having a value less 
than zero and divide by J . Moreover, the empirical cumulative distribution 
function (CDF) of Z can be plotted on probability paper. 

It is a well known fact that the Monte Carlo solutions tend to be 
expensive relative to the accuracy obtained. Often a large sample size 
is required, particularly in structural reliability analysis when the 
values of the probability of failure are typically small. For example, 

A 

it might require 10 samples to provide a reasonable estimate of probability 

-3 

of failure in the range of p^ > 10 . If the failure function is compli- 
cated so that a costly analysis is required for each Z^, Monte Carlo be- 
comes impractical. Nevertheless, the method is easy to apply, able to 
solve complicated problems, and therefore is a valuable tool for design 
and for research. 

A limited Monte Carlo analysis is employed later to check the results 
of the proposed method. 

FAST PROBABILITY INTEGRATION TECHNIQUE 

The Rackwltz-Fiessler (R-F) algorithm, used herein, is well documented 
in literature [9, 12, 13, 14]. The lollowlng procedures are the summary 
of the method. 

1. Transformation of the basic variables to the reduced variables, 
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where y i end o, are the jean and the atandard deviation of 
respectively. Note that u^ are random variables having a mean 
of sero and atandard deviation of unity. 

2. Determine the minima distance! 6. froa the origin to the limit 
state on the reduced coordinates. 

3. Obtain the "equivalent noraal" values of and o^ for those 
random design factors which do not have a normal or lognormal 
distribution. This Is done by Imposing the conditions that the 
density and distribution functions of the equivalent normal 
are the same as the basic variable at the design point. 

4. Repeat steps 1. to 3. until 0 converges. 

The point on the limit state closest to the origin is called the design 
point 0 is called the safety index and is considered to be a measure of 
reliability because It is inversely proportional to the probability of 
failure. 


It is a common practice to calculate 0 by Iteratively solving the 
following equations [9, 14] 



u’ - 0X i (13) 

g(x 1 e, x 2 e, . . . x n 0) - o (14) 

where the partial deri\ates are calculated at the design point, jj'j 
are the direction cosines of the vector that define the design point. 
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Note that this scheme requires computation of Vartial derlvatlvea of the 
failure function, g. In the case when the failure function cannot be ex- 
pressed explicitly, the method becomes inefficient. Furthermore, in many 
practical situations the failure function may not be explicitly differenti- 
able, and it is necessary to use numerical dlfferention. Again the failure 
function needs to be known. 

B can also be computed by solving the constrained optimatlon problem: 



(15) 


gty) * 0 


(16) 


The computations of partial derivatives are avoided (e.g., see Ref. 11). 
For a complicated failure function, it is suggested that this scheme be 
adopted and the failure function be approximated by a closed form alge- 
braic expression so that 8 and j^’can be effectively obtained. 

A probability of failure computed by 


" t(-B) 


(17) 


provides a reasonable estimate of the "exact" value of p^ for many cases 
[A, 16]. A study of the R-F algorithm has demonstrated that for many 
cases, (a) estimate of p^ are "good," to "excellent" and (b) R-F computer 
costs are a fraction of Monte Carlo, e.g., 1/50 for some examples fl6]. 


POLYNOMIAL APPROXIMATION TO THE FAILURE FUNCTION 

If the failure function can be approximated by a simple (e.g. , poly- 
nomial) equation, B can be determined effectively by solving the optimiza- 
tion problem (F-.s 15 and 16). The entire process of computing 3 and then 
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p . is summarized by the block diagram of Fig. 1. The method will be de- 
veloped and demonstrated by solving a sample problem. 

Assume, in general, that there are n + 1 variables involved in a 
failure function expr sion 

g(U , U, , . . . U ) - 0 (18) 

o l n 

Assume also that the failure function is differentiable with continuous 
partial derivatives in the whole domain. By expanding it at a point 
a = (a , a .... a ) by Taylor series, the following polynomial results: 

1 O 1 T ' 


*<%) 


sty) + 1 

i=0 







a i ) 


2 


n-1 n 
+ 22 I 

i-0 j=i+l 


» 2 x 


3U i 3U J 


(U. - .,) (Uj 

A 


- a^) + H.O.T. (higher order terms) 

(19) 


In order to compute 8 using the optimization scheme of the R-F algorithm, 

it is necessary to rewrite the equation such that any one of the variables 

is a function of the remaining variables. Thus Eq. 18 becomes, when one 

chooses U , 
o 

U - h (U. , U_, . . . U ) (20) 

o X I n 


such that the failure function is 

gtl P - D 0 - h (0 1§ U 2 . . . . U n ) - 0 (21) 

By substituting Eq. 21 into Eq. 19, it follows that 
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FIGURE 1. Block Diagram for Computation of the Safety Index , B, when 
the Failure Function is Not Available as a Closed Algebraic 
Form. 


Specify the random design factors. 

(The distribution types and parameters.) 


Select several points for each variable 
where a solution is desired. The points 
should embrace the design point, but one 
has to guess where that will be. 


Obtain a solution, using computer analy- 
sis , at each one of the selected point 
sets. 


Curve fitting: Use a polynomial to 

approximate the failure function. 
The limit state is thereby defined. 


Improve the accuracy of S by 
selecting more appropriate 
points using design point values. 


The R-F algorithm is applied to com- 
pute the design point and safety 
in«ex, B. 


The estimate of the probability of 
failure is = 4>(-B) 
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* 


n-1 n 
+ 21 E 
i-1 j-i+1 


a 2 h 




<o i - “i >< n 3 


a^) + H.O.T. 


( 22 ) 


Simplifying, 


U 

o 


n 
+ I 
i-1 


(«!- 


a l> 


n 
+ I 
i-1 


n+i 


(»i- 


».) 


n-1 n 
+ E E 
i=l j-i+1 


a (n-|)(i+l)+j 



- a t )(U 


j 


- 3^) + H.O.T. 


(23) 


where the ct's are defined by comparing Eq. 22 and Eq. 23, Eq. 23 can be 
further simplified as 


k 

U - a + E a x. + Error (H.O.T.) (24) 

o o ^ . i i 

where 

^ - a^ for i - 1, n 

-(U - a^) 2 for i - n + 1, 2n 

= (U. - a. )(U - a ) for i - 2n + 1, k - t(n + 2)(n + l)/2] -1 

1 1 J J 

This is the familiar multiple linear regression model where are inde- 
pendent random variables. (Note: The model is linear in the parameters 

(the a's) regardless of the shape of the surface it generates.) The least 
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square normal equation may be written in matrix notation aa (see e.g. 

Ref. 10) 

[X) T [X] [a] - [X] T [U o ] (25) 

where 

(26) 


[a] - 

Note that there are k + 1 normal equations, one for each of the unknown 

a., m is the number of observations or data sets. Given the data sets x . . 
i lj 

and (U )^, the solution will be the least square estimators, a^. 

SELECTION OF THE DATA POINTS 

Assume that the failure function is to be approximated by a truncated 
polynomial equation. In a general data analysis problem, data obtained 
from experiments are subject to random errors. A statistical analysis 

A 

of estimators are required to validate the model. In such cases, 
accuracy of the model will be heavily dependent upon the sample size, m. 


a 

o 

°1 

, [U ] = 

<»o>l ' 

«o>2 

• 

• 

o 

• 

• 

. v 


1 

0 

• O 

J=> 

1 


IX) 


1 *11 x 12 • • • x lk 


. x 21 x 22 . . . x 2k 


• • 


1 X , X „ . . . X , 

ml m2 mk 
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In this study, however, the data sets x ^»(U q )^ are those values which 
satisfy the exact failure function. Therefore, smooth functions are 
anticipated and statistical analysis is not necessary. The number of 
data sets for curve fitting can be chosen relatively small. More pre- 
cisely, the minimum number of data sets (m) is equal to the number of 

V 

The key to constructing an effective approximating model is to choose, 
for curve fitting, data points (centered around the expanding point 
from a domain that is relatively small to give satisfactory agreement 
between the exact and the approximating limit state, but wide enough to 
cover the unknown design point. When R-F algorithm is used, the model 
need be accurate only in the region of the design point rather than the 
whole domain. If the expanding point, which must be chosen at the 
outset, is relatively close to the design point, "accurate" results can 
be expected. 

In the general case, the design point may be difficult to predict in 
advance. However, it is reasonable to assume that design point values of 
stress variables would be greater (may be much greater) than its mean or 
median. Similarly, design point values of strength variables would be 
less than its mean. Hence, one can at least make a reasonable first guess 
of the neighborhood of the design point. Another observation which is 
useful is that the failure functions are expected to be relatively smooth 
and unlikely to have discontinuities or other significant nonlinearities 
in the neighborhood of the design point. Therefore, it seems reasonable 
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to use a quadratic polynomial model to approximate the failure function. On 
the other hand, It is worth mentioning that the simple relationship be- 
tween p f end 6 of Bq. 17 is a good approximation when the limit state sur- 
face can be approximated by a tangent hyperplane at the design point 
[5]. In such cases it Is evident that a quadratic model will be suffici- 
ent to approximate the limit state in the region close to the design point. 

As an example of solving using minimum of data points, consider 
the case where there are four design variables involved. If the quadratic 
polynomial without Interaction terms (or mixed terms, e.g., (U^ - a., ) . 

(U 2 - a ? ) , . . . etc.) is chosen, the model is 


U o ■ “o + V“l • V + “2 (U 2 - a 2> + °3 (B 3 * a 3> 

+ a 4 (U 1 - Oj) 2 + a 5 <U 2 - a 2 ) 2 + 0 6 <Hj - .j) 2 (27) 


The expanding point ,a has 

been 

chosen 

and it is required 

to 

compute the 

estimates, a . . . a, . 

o 6 

These 

estimates can be computed 

by 

calculations 

relating the design parameters 

at the 

following specific 

points. 

h 



S 

U 

o 



a l 

a 2 


a 3 

<Vi 



a ± + AUj 

a 2 


a 3 

(U o)2 



3 1 - AU 1 

a 2 


a 3 

(U o>3 



a l 

a 2 + 

4U 2 

a 3 

<“o>4 



a l 

a 2 “ 

iU 2 

a 3 

<"o>5 



a l 

a 2 


a 3 + AU 3 

(U o>6 



3 1 

a 2 


a 3 ' AU 3 

< U o>7 
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The AU^'s must be chosen. They must be large enough to embrace the design 
point, but not too large so as to Introduce significant errors In the poly- 
nomial approximation. Here the analyst must again use his Intuition. 

Finally, the estimates of a's are computed by noting the following. 

^ A A 

ct Q is the value of U q at the expanding point; ct^ and can be determined 

A 

when U 2 - a 2 and U 3 ■ a 3 , . . . etc. Thus, can be determined from the 
above seven data sets by choosing appropriate values of ± AU, which define 
the range of data for each variable. 

A 

In general, one can solve using matrix operations (Eq. 25), where 
x matrix for the example is 


[X] = 


1 

1 

1 

1 

1 

1 

1 


0 0 

AU X 0 

-Al^ 0 

o au 2 

0 _AU 2 
0 0 

0 0 


0 0 

0 AuJ 

0 AU2 

0 0 

0 0 

au 3 0 

-AU 3 0 


0 0 

0 0 

0 0 
AU* 0 
AUj 0 
0 AU* 
0 AU2 


(28) 


If the interaction terms are included 

A A 

sets are required to determine to a^. 


, a minimum of three more data 
As an example. 


U 1 


a 

a 


1 

1 


+ AU, 


- AU, 


0 



U 3 

u 

0 

- AU 2 

0 

(U o>8 

0 

a 3 + “3 

o\ 

0 

+ AU 2 

a 3 - “3 

0 

H 

O 

£3 


Similar tables can be easily constructed when the number of design variables 
(n + 1) is greater than four. In general, the minimum number of data sets 
required is 2n + 1 if the interaction terms are excluded, and is 
(n + 2)(n + l)/2 for a complete quadratic equation. 
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When using the minimum number of date points to generate the model, 
the values at the data points of the model are exact. Then the curve 
fitting procedure described above is actually an effort of trying to 
approximate the exact surface by interpolation. To determine suitable 
values of the three data points for each variable, namely, the expanding 
point and the other two points which define the "interpolation range," 
consider the following example. Assume that a stress random variable 
U has a known distribution, but no information is available about the 
value of the design point. As shown in Fig. 2, a decision is made to 
choose the Interpolation range from the median of U (data point 1) to the 
value of U where its CDF is 0.9999 (data point 3), and let the expanding 
point to be in the middle (data point 2). For a strength variable, points 
2 and 3 should be chosen at the left tail region. In doing this, it is 
anticipated (as a first approximation) that the design point will lie with- 
in the interpolation range. To improve the approximation, a new range 
may be chosen if the design point lies outside this range in the first 
reliability computation by the fast probability integration technique. 

As demonstrated later, it is thought that one or two tries will be suf- 
ficient for generating satisfactory results. 

In the first attempt, in general, one would like to choose a wide 
Interpolation range to ensure that the range will cover the design point. 
However, it is obvious that a wider interpolation range may result in 
poorer approximation. In addition, if significant nonlinearities exist 
in the range, the quadratic curve or surface may turn into an inappropri- 
ate direction very fast which may result in a problem of obtaining 
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FIGURE 2. An Example of Choosing Data Points for Generating Curve Fitting 
Model (Stress Variable Case). 
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unreasonable design point. It Is helpful to have a good physical feel 
of the problem. Consider the example of fatigue reliability analyala 
mentioned before. For e given life N, fatigue strength (defined by e) 
should increase if cl being Increased. In other words, 3e/3o| should 
always be positive. Suppose a linear polynomial model Is chosen. As a 
first approximation, the data points (e^, o£, . . . ) will be fitted by a 
straight line on the two dimensional coordinates plane (e.g. , e - 
plane). But the data ere exact, and It follows that the sign of 3e/3o^ 
will be preserved. Hence, the linear model Is excellent In preserving 
the basic functional relationships between strength and stress variables. 

In general, however, the quadratic model is more appropriate because 
of the wide range of good fit of the failure surface. But It is also 
clear that the slopes or (first) partial derivatives of a second degree 
equation can be positive, negative, or zero. Therefore, when a "strictly 
monotone" nonlinear function is approximated by a quadratic equation, It 
is no longer a strictly monotone function. As a result, there exist 
certain domains in which a strength variable will behave like a stress 
variable, and vice versa. Consequently, there is a possibility that the 
design point will converge to the undesirable domain where very poor ap- 
proximation of the limit state is unavoidable, resulting in unreasonable 
design point. Therefore, it is important to examine the model before 
calculating 6. 

Instead of checking the result of the curve fitting by plotting, it 
is more efficient to calculate the "effective domain" of the model where 
the basic functional relationships mentioned above are preserved. As 
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Illustrated in Fig. 3 the effective range can be date rained for each 
variable by calculating the value at which > partial derivative la aero. 
For instance, recall Eq. 27 for a quadratic polynomial with no Interaction 
terns. Consider U o as the dependent variable, and U^, Uj, as the in- 
dependent variables. Fron the conditions of 3 U q / 3U^ - 0, the "critical 
points" of denoted by U*, are 


U, 


a l 

r2a 4 

(29) 

a 2 

2 2<* 5 

(30) 

a„ 


3 

3" 2a 6 

(31) 


Thus, the effective range can be easily determined. Similar reaults can 
be immediately Inferred when more design variables are involved. 

By Introducing interaction terms to form a complete quadratic equa- 
tion, Eq. 29 is replaced by 

U 1 “ (a l " *T * _ 2^7 (U 2 " a 2* " 2^7 (U 3 “ *3 ) (32) 

4 4 4 


Considering the extremes of the interpolation range where U 


*i * 4U i 


for and U^, Eq. 32 becomes. 


(a i - 5S7 > * \ « 4 V - i (t 

4 4 4 


(33) 


therefore, at least at one of the extremes, will be 


u t ' (a i * S7 > - 

4 


a 7 ... 


a 

CD 

2sr iU 2 

4 


IT 4U 3 
4 


(34) 


which is at best, equal to the value obtained from Eq. 29. Similar argu 

* 

monts can be used for the other U^, and for the case of more variables. 


22 



FIGURE 3. Definition of the Effective Range for the Quadratic Polynomial 
Curve Fitting Model. 
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Thus, while the Introduced Interaction terns have the proarimes of inprov- 
ing the accuracy of the model, the effective range becomes smaller, 

A general rule regarding the quality of the polynomial expression 

* 

is that the model la good if all located well outside the interpola- 

tion range (See Fig. 3) . It is clear that the effective range will be 

much wider than the Interpolation range if the failure function is approxl- 

* 

mately linear in the interpolation range. On the other hand, when lies 
Inside the data range, it is implied that significant nonlinearity exists 
in the interpolation range. A "bad fit" results in some region of the 
interpolation range. Likely the fit is worse In the "extrapolation range." 

Using the information of effective range, there are several ways of 
overcoming this "nonlinearity problem": (a) Restrict the value of the 
variable within its effective range. This will be demonstrated in the 
following example of the fatigue life reliability analysis; (b) Use linear 
model as a first approximation followed by a quadratic model using more 
appropriate Interpolation range; (c) Use a piecewise (splines) model by 
splitting the interpolation range. In the limiting case, any curve can be 
considered as connected by a series of straight lines. 

The methods atove are useful when a range of p^ need be estimated, 
and a single quadratic curve is unable to provide a good fit in the 
relatively wide domain because of the nonlinearity of the failure function. 
For the case of single design-point checking problem, "range reduction" 
can be very effective simply because a narrower interpolation range always 
provides a better effective range. 

In conclusion, it is thought that generally the quadratic polynomial 
curve fitting is an appropriate method to approximate the failure function. 
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Furthermore, because the interaction terms tend to reduce the effective 
range and more data point sets are required for modeling, it is recom- 
mended that a quadratic polynomial model without the Interaction terms 
may be adopted for an efficient estimation of p^ and the design point 
values. Using a minimum amount of data points, this model will be 

applied in the following example to demonstrate its effectiveness. It 

* 

will be shown that critical values, are useful for validating and 
improving the model. Finally, it is believed that accurate 0 can be 
effectively obtained through the "design-point searching" and the 
"interpolation-range reduction" processes regardless of the complexity and 
nonlinearity of the failure function. 


EXAMPLE - FATIGUE LIFE RELIABILITY ANALYSIS 

Demonstrated here is v.he general method of producing a quadratic 
polynomial model where interaction terms are neglected. The sample problem 
involves fatigue reliability using local strain analysis to predict life. 

A zero-to-peak constant amplitude load producing the nominal stress history 
as shown in Fig. A is applied to the notched member. The peak stress S 
is not known exactly, and is considered to be a random variable. An il- 
lustration of the basic equations of local strain analysis (Eqs. 5 to 7) 
are provided in Figs. 5 and 6. Fig. 5 shows that a and e are determined 
from Neuber's rule and the cyclic stress-strain curve. Fig. 6 illustrates 

the relationship between strain amplitude at the notch (e ) and cycles 

a 

to failure (N) for a given and r.^. Both and are considered to be 
random variables, describing the scatter in fatigue data. 

Three different applied stress conditions are examined using the pro- 
posed method described herein. The results are compared with those using 
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FIGURE 4. Nominal Stress Applied to Notched Member. 
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Strain Amplitude, 



FIGURE 6. Strain Life Curve for Case A (Using Median Value of S, o' { and e^). 
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Monte Carlo method where 100 samples are generated In each case to plot 
the empirical CDF. 


Case A : 

The material analyzed is Uaspaloy B at 1000° F. The median nominal 

'V, 

stress is S * 310 MPa (45 ksi) , The data is tabulated in Table 1. To 
obtain wide range of effectiveness of the polynomial fit, data points 
are chosen as described in Fig. 7. 

Cycles to failure N is computed using a local strain analysis program. 
Seven data sets for curve fitting (analysis is repeated seven times) are 
given in Table 2. The scale of N has been changed for easier data handling. 
Once the data sets are generated, any variable can be treated as the de- 
pendent variable. In this example, S was chosen. 


S * “o + “l <Cl f " a l ) + “2 (e f • V + “3 (log 8 - °3> 

+ (oj - a^ 2 + Oj (cj - a 2 ) 2 + (log N - a 3 ) 2 

where the expanding point values, a^, can be found in Table 2. By the 
transformation of ” a i» e *c., as defined in Eq. 24, S can also be 

written as 


S ■ a + a, x, + ... a, x, 
oil 66 


(36) 


Therefore, can be solved simultaneously from Eq. 25 where the X matrix 
can be constructed using data in Table 2. The result is: 


S = 497.7 + 0.268 (o' - a ) + 20.8 (e’ - a 2 ) - 168 (log N - a 3 ) 

+ 0.000027 (o' - a 1 ) 2 - 4.52 (-J - a 2 ) 2 + 37.5 (log N - a 3 ) 2 


(37, 
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TABLE 1. 
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Data Information for Example 


(a) Random Variables 


Parameter 

Distribution Type 

Median 

Coefficient of Variance 

S 

lognormal 

310.2 MPa* 
(45 ksi) 

0.20 

«£ 

lognormal 

1841 MPa 
(267 ksi) 

0.05 


lognormal 

3.47 

0.43 


S - 414 MPa (60 ksi) for Case B, $ « 551 MPa (80 ksl) for Case C 


(b) Constants 




Parameter 

Value 

Parameter 

Value 

K t 

3.0 

s 

0.108 

b 

-0.0843 

c 

-0.9126 

E 

206820 MPa 

A 

1518 MPa 


(30,000 ksi) 


(220.4 ksi) 

Note: o - 
o 

0 



■ - .. . . --I 
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TABLE 2. Data Point Sets for Curve Fitting 
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Set Number 

S, in MPa 

oj.in MPa 


i° 8 io N * in c y cles<b) 

l a 

497.7 

1 

2.069 

3.859 

2 

310.2 


2.069 

6.191 

3 

685.1 


2.069 

2.935 

4 

497.7 

1507 

2.069 

3.611 

5 

497.7 

1841 

2.069 

4.149 

6 

497.7 

1674 

0.668 

3.644 

7 

497.7 

1674 

3.470 

3.984 

Interpolation 

Range 

310.2-685.1 

1507-1841 

0.668-3.47 



2.935-6.191 

0 

Expanding point 

^Computed using a local strain analysis program originally written 
by W. R. Brose. 



32 















Now examine the quality of the polynomial expression. Consider the 
effective range for each Independent variable. Using Eqs. 29 to 31, the 
values at 3S/3o* - 0, etc., are 

(oj, e£, log N)* - (-3290, 4.4, 6.1) 

The value of (op* Is negative implying the nonlinear effect of is 
very small. The value of (cp* when compared with the interpolation range 
of e| shown in Table 2 also suggests that the nonlinear effect of e£ is 
not significant. Thus, the quadratic model seems to be satisfactory for 
and in the chosen ranges of interpolation. The value of (log N)* 
is close to the interpolation range implying a moderate nonlinear effect 
of log N in its chosen range. 

The R-F algorithm was used to compute p^ for several values of N^. 
The results are shown in Table 3. It is noted that the highest value of 
log N q is chosen to be 5.5, a value less than (log N)*, 6.1. This point 
is inside the Interpolation range and reasonable p^ estimates are expected. 
There is no limitation for choosing lower log N q . The values listed in 
Table 3 are in the range where comparisons can be made with Monte Carlo 
results, i.e., where the lowest possible estimate of p^ is about 12. 
Observing that all the estimated design point values in Table 3 are well 
within the interpolation range shown in Table 2 one notes that perhaps 
accuracy could be further improved by reducing the interpolation ranges. 

Monte Carlo was used as an independent check on the quality of results 
of the proposed method. Unfortunately, each trial requires running the 
local strain analysis program with a random sample of S, a ’ , and ep This 
process is extremely costly and inefficient. Therefore, only a sample of 
size 100 was obtained. 
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TABLE 3. Results of the Reliability Analysis for Case A Using the 
Proposed Method 



Design Point 

Probability of 
failure, X 

N q (cycles) 

l°gN o 

S, in MPa 


E f 

3,160 

3.5 

605 

1779 

2.88 

0.026 

10,000 

4.0 

520 

1785 

2.95 

0.34 

31,600 

4.5 

454 

1793 

3.04 

2.2 

100,000 


407 

1802 

3.13 

7.3 

316,000 


378 

1810 

3.21 

14.2 

- - ■ 
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Computed values of from the R-F algorithm ere plotted In Fig. 8. 

The empirical CDF generated from Monte Carlo method is also plotted on the 
same lognormal probability paper. Assuming that the p^ of the present method 
are accurate, the upper and lower 1Z points for the maximum likelihood 
estimator (MLE), p, of p are given. The calculation is similar to that 
of a confidence Interval for the binomial. If the circled points in Fig. 8 
are accurate, then there is a 98Z chance that p obtained from an experiment 
would lie in the interval shown. Even though the algorithm used to esti- 
mate Pj was really (i - y)/n, not exactly the MLE, the evidence of Fig. 8 
suggests good agreement between the present method snd Monte Carlo. 

The median of S in this case was sufficiently small so that generally 
the deformation of the material is In the elastic range. However, because 
S has a large coefficient of variance, some plasticity is expected for 
samples having a higher than average stress. The platidty effect (c^ term) 
could be Important in highly stress specimens. Hence, the lower tall of the 
CDF of N is not a straight line as it would be on lognormal paper if S 
were deterministic. This exercise demonstrates that extrapolating p^ using 
Monte Carlo results will likely produce Inaccurate estimates of p^. This 
is a principal reason why the number of simulations should be larger when 
the Pf to be estimated is smaller and why Monte Carlo can be so expensive. 

Case B 

The data is the same as in Case A (See Table 1) but the median of 
nominal stress is Increased to 414 MPa (60 ksi) so that more cyclic plas- 
ticity is expected at the notch. 
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FIGURE 8. Reliability Analysis of Fatigue Life (Case A); A Comparison of the 
Present Method with Monte Carlo. 
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The data points for S ara (414, 663, $13) MPa. The aethod for com- 
puting these was described in Fig. 7. Following the same procedures as 
in Case A, seven data points were obtained. The expanding point of log N 
is 3.0, and the data of log N ranges froa 2.5 to 4.7. 

The result of the curve fit, choosing S as the dependent variable, la 

S - 663 + 0.217 (o! - a.) + 106 (ei - a-) - 425 (log N - a,) 

(38) 

+0.000074 (o’ - a x ) 2 - 34.1 (ej - a 2 > 2 + 164 (log H - a 3 > 2 
The values at 3S/3 - 0, etc., are 

(o* . ej, log N)* - (207, 3.6, 4.3) 

Again, both (op* and (c£)* values are well-outside the interpolation 
range shown In Table 2 implying satisfactory effective ranges for o’ and c^. 

However, (log N)* is Inside the Interpolation range (2.5, 4.7), 
suggesting the possibility of a poor fit in some regions because of sig- 
nificant nonlinearity between S and log N in the interpolation region. 

Therefore, in forming the limit states of Eq. 38, values of log N q should 
be limited to, say 4.0 as higher errors uay result for larger log N q . 

The design points and p for several values of log were computed 
using the same method as described for Case A. The results of this first 1 

approximation is shown in Table 4(a). Values of the Pg are plotted on j 

Fig. 10 as the dashed line. j 

The fact that the (log N)* value lies inside the interpolation range j 

implies that the degree of nonlinearity between S and log N is significant 

i 

i 

in this range. A higher degree polynomial equation is required to 1m- i 

i 

i 

prove the fit unless the range can be reduced. Suppose, for instance, 

the goal of the reliability analysis of log N ranges froa 3.0 to 4.5. 

o 
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The design point information generated from the first approximation 
(Table 4(a)) reveals that the interpolation ranges may be reduced to gen- 
erate improved, second approximation. It should be emphasised that when 
only point estimate is required (e.g., log - 3.0), the interpolation 
range for the second approximation can be made very small so that the 
quadratic model becomes a good approximation. In this study, however, 
the narrower interpolation range is chosen to be wide enough to cover the 
region of log N q of Interest in this study. 

Using the design point Information in Table 4(a), the new and reduced 
interpolation ranges are determined as shown in Table 4(b). The improved 
quadratic model becomes: 

S « 551 + 0.212 (ol - a.) + 19.3 (ei - a,) - 165 (log N - a,) 

t l t i i (39) 

+ 0.000074 (o f ' - a T ) 2 - 4.08 (ej - a 2 ) 2 + 55.1 (log N - a J ) 2 

where the expanding point values can be found in Table 4(b). The values 
at ;>S/r»o* ■ 0, etc., are 

(o', e', log N)* - (390, 4.6, 5.2) 

By comparing these values with the interpolation ranges shown in Table 
4(h) It is found that the effective ranges for al? three variables are 
satisfactory, suggesting that the model has been improved. To demon- 
strate the Improvement made in the model, the exact and the approximate 
S - N curves are drawn in Fig. 9 for both the first and the second approxi- 
mation (Note: Values of o ^ and c ^ are at the expanding point). The figure 

clearly shows that the model is much better when (log N)* lies outside the 
interpolation range suggesting that (log N)*, etc., which defines the ef- 
fective range, is a convenient index of "measuring" the quality of the 
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TABLE 4. Results of the Reliability Analysis Using the First Approximating 
Model, and the Data Range for Generating Second Approximating 
Model (Case B). 


(a) Results Using First Approximating Model 




Design Point 

Probability of 

N q ( cycles) 

iogN, 

S, in MPa 



failure, X 

1,000 

3.0 

658 


1.75 

0.19 

1,778 

3.25 

587 . 


1.95 

1.14 

3.162 

3.5 

539 

... . .■ 

2.20 

3.97 

5,623 

3.75 

508 


2.47 

8.93 

10,000 

4.0 

490 

1823 

2.69 

14.2 



(b) Data Range for Generating Second Approximating Model 


L ’ 1 

Variables 

Data Range 

Expanding Point 

S, in MPa 

655 - 448 

551 

0 ^, in MPa 

1806 - 1834 

1820 

*i 

1.75 - 2.80 

2.275 

a log N 

3.15 - 4.58 

3.686 

a Corresponding Values 



39 



















Applied Stress 



log N 


FIGURE 9. Comparison of the Exact and the Approximating Failure Function. 



model, and the degree of nonlinearity In the failure function. 

The new design points and the associated using Eq. 39 are shown in 
Table 5. The design point values of as well as the value of S when 
log N q - 3.0 are outside the data ranges indicating the extrapolating values 
of Eq. 39 are being used. However, these values are well within the ef- 
fective ranges and closer to the interpolation range. Therefore, it is not 
unreasonable to assume that good approximation can still be obtained. 

Values of p^ of the improved model are also plotted on Fig. 10. 

Following the same procedure as in Case A, a 98Z interval on the estimates 
of p^ is constructed which covers the Monte Carlo results. The improve- 
ment of the model seems justified. By using more appropriate data ranges 
or splitting the data range into two parts, accurate Rackwitz-Fiessler values 
of 8 can be obtained for a wide range of fatigue life. For a specified 
single value of design life, it is believed that no more than two tries 
are necessary to obtain satisfactory results. 

Case C 

For the previous examples, the median nominal stress is further increased 
to 551 MPa (80 ksi) so that even more cyclic plasticity is expected. 

Otherwise, the parameters are the same as in Cases A and B. 

Applying the same methods as described in previous cases, the data 
points for S are (551, 885, 1218) MPa. The corresponding expanding point 
value of log N is 2.55, and the interpolation range of log N is from 2.05 
to 3.48. The result of the curve fit is: 
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TABLE 5. Improved Results of the Rellabilit/ Ant* lysis for Case B Using 
the Proposed Method. 



Design roint 

Improved 

First 

T? a a 4 a m «* 

N q (cycles) 

!og N q 

S, in MPa 



P f , * 

LalJlialcS 

of P f , x 

1,000 

3.0 

700 

1808 

3.05 

0.34 

0.19 

1,778 

3.25 

645 

.810 

3.08 

1.12 

1.14 

3,162 

3.5 

596 

1814 

3.13 

2.96 

3.97 

5,623 

3.75 

555 

1817 

3.18 

6.43 

8.93 

10,000 

4.0 

520 

1821 

3.23 

11.7 

14.2 


4.25 

493 

1825 

3.28 

18.2 

— 


4.5 

472 

1828 

3.32 

24.7 

— 
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FIGURE 10. Reliability Analysis of Fatigue Life (Case B) ; A Comparison of the 
Present Method with Monte Carlo. 
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S - 885 + 0.133 (Oj - a x ) + 210 (e* - * 2 ) + 702 (log N - « 3 ) 

+ 0.000012 (aj - a^ 2 - 77.1 (ej - a 2 > + 370 (log N - a 3 > 2 


(40) 


The values at 3S/So^ - 0, etc., are 

(o£, ej, log N)* - (-3870, 3.43, 3.5) 

Thus, satisfactory effective ranges for and e’ result. The value of 
(log N)* Is approximately equal to the value at right end of the interpola- 
tion range suggesting that the range of interpolation of S is too wide. 

To improve the model, a decision can be made to reduce the interpolation range 
either before or after the design points are estimated using Eq. 40. The 
method of improving the accuracy through range reduction process has been 
demonstrated in Case B, therefore it is not repeated here. 

To obtain the first estimation of reliability analysis of N, several 
values of log N q are chosen to form the limit states. The highest value is 
restricted to 3.25 which is somewhat less than the value of (log N)*. The 
design points and p^ using R-F algorithm are shown in Table 6. The values 
of e’ and indicate that cyclic plasticity indeed plays the dominant role 
in the region of high p^. p^ values are plotted in Fig. 11. A 98% interval 
on the estimates of p^, and the Monte Carlo results are also plotted. It 
is noted that for this high stress case, lower tail CDF of should approach 
a straight line on the lognormal paper if the elasticity effect is being 
neglected. Thus, the result of the first estimation is reasonable. 


SUMMARY 

A numerical procedure is proposed for estimating risk for the case 
when the failure function is not available as a closed algebraic form. 
The safety index is computed by the R-F algorithm after the 
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TABLE 6. Results of the Reliability Analysis for Case C Using the Proposed 
Method. 




Design Point 

Probability of 

(cycles) 

log N o 

S, in MPa 


e f 

failure , X 

178 

2.25 

927 

1825 

1.25 

0.015 

316 

2.5 

786 

1828 

1.46 

0.28 

562 

2.75 

706 

1831 

1.75 

1.8 

1000 

3.0 

663 

1833 

2.06 

5.7 

1778 

3.25 

641 

1834 

2.30 

10.5 



\ 
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failure function la approximated as a simple polynomial equation. The proposed 
method is described by the following step by step procedures (refer to Fig. 1): 

1. Specify the random design factors. 

2. Make a somewhat conservative but reasonable first guess of the design 
point to define the range of Interpolation. Select several points 
from the range for each variable where a solution Is required. 

3. Solve the failure function, using computer analysis, to obtain a 
solution at each one of the point sets selected. The data sets 
for curve fitting are thereby provided. 

4. Obtain a polynomial equation using a curve fitting technique. It 
is recommended that a quadratic equation without the Interaction 
terms may be adopted for an efficient estimation. 

5. Check the quality of the polynomial fit (e.g., define effective 
ranges) as an aid in forming the limit state where a reasonable 
agreement between the exact and the approximating equation can be 
anticipated. A reduction of the ranges of the selected points 
may be required if the failure function is highly nonlinear in 
the selected interpolation range. 

6. Apply R-F algorithm to the approximating limit state to compute 
the design point and 8. 

7. Compare the design point values with the ranges of the selected 
points. To improve the approximation in the region of design point, 
go to Step 2, but adjust the values and ranges of the points for 
each variable. In general, accuracy can be improved by reducing 
the ranges which should still embrace the design point. 

8. Estimate the probability of failure using p f - d> (—3) . 
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